新的ngs流程该如何学习之m6A学习大纲
上个月生信入门的学员花了一周学习了m6A的生信相关知识,汇报给了我一个笔记,主要记录m6A的建库原理和分析流程。跟我前些天提到的教程:新的ngs流程该如何学习(以CUT&Tag 数据处理为例子),有异曲同工之妙,所以分享给大家。确实我自己是不太可能把所有的ngs流程全部录制视频的,但是我仍然是希望传达学习方法给到大家。
如果你看过我表观组学,比如《ChIP-seq数据分析》 和 《ATAC-seq数据分析》 就会发现其实这个m6A数据处理大同小异的,当然了,肯定是会有一些细微差异是需要注意的。
首先可以参考其它NGS组学技能
比如我在生信技能树就录制了多个NGS组学数据分析教学视频的,都免费在B站,而且组建好交流群,见:
免费视频课程《RNA-seq数据分析》 免费视频课程《WES数据分析》 免费视频课程《ChIP-seq数据分析》 免费视频课程《ATAC-seq数据分析》 免费视频课程《TCGA数据库分析实战》 免费视频课程《甲基化芯片数据分析》 免费视频课程《影像组学教学》 免费视频课程《LncRNA-seq数据》 免费视频课程《GEO数据挖掘》 肿瘤基因测序
如果你认真看了我是五年前是如何设置这些ngs教学大纲就应该是很容易明白,一个新的ngs流程该如何学习。无非是找案例,看综述,画流程,安装软件,调试测试数据,解释图表,最后延伸到自己的数据。
但是,什么是一个入门级别工程师呢,我认为首先必须具备生物信息学编程基础知识,就是基于R语言的统计可视化,以及基于Linux的NGS数据处理。
必须生物信息学编程基础知识
再怎么强调生物信息学数据分析学习过程的计算机基础知识的打磨都不为过,我把它粗略的分成基于R语言的统计可视化,以及基于Linux的NGS数据处理:
把R的知识点路线图搞定,如下:
了解常量和变量概念 加减乘除等运算(计算器) 多种数据类型(数值,字符,逻辑,因子) 多种数据结构(向量,矩阵,数组,数据框,列表) 文件读取和写出 简单统计可视化 无限量函数学习
Linux的6个阶段也跨越过去 ,一般来说,每个阶段都需要至少一天以上的学习:
第1阶段:把linux系统玩得跟Windows或者MacOS那样的桌面操作系统一样顺畅,主要目的就是去可视化,熟悉黑白命令行界面,可以仅仅以键盘交互模式完成常规文件夹及文件管理工作。 第2阶段:做到文本文件的表格化处理,类似于以键盘交互模式完成Excel表格的排序、计数、筛选、去冗余,查找,切割,替换,合并,补齐,熟练掌握awk,sed,grep这文本处理的三驾马车。 第3阶段:元字符,通配符及shell中的各种扩展,从此linux操作不再神秘! 第4阶段:高级目录管理:软硬链接,绝对路径和相对路径,环境变量。 第5阶段:任务提交及批处理,脚本编写解放你的双手。 第6阶段:软件安装及conda管理,让linux系统实用性放飞自我。
m6A建库原理及分析流程学习
实验环节获取RNA
提取RNA,注意避免RNA的降解,RIN>=7.5 , 从左至右三个峰分别是5S、18S、28SrRNA,上面的图有很多杂乱小峰,说明有RNA降解,质量不太好;下面的图RNA质量较好.
建库测序
RNA打断片段化(如何片段化?后续查找学习)
IP文库:m6A 抗体富集m6A RNA
Input文库(相当于转录组文库)
m6A修饰区域有Reads的堆叠.
分析
数据质量可视化
原始下机数据
FASTQ
1.Read ID
2.Read序列
3.一般为一个"+"
4.碱基质量值
Q=-10lg perror 值越高质量越高
可视化软件:fastqc/multiqc/fastp
测序前面几个碱基不均衡,需去除;因为m6A含量低,所以接头含量高
数据过滤
软件:TrimGalore/Trimmomatic(https://www.jianshu.com/p/7a3de6b8e503)
Input:Raw FASTQ
Output:Clean FASTQ
质量比较好的fastq过滤后的数据Q20>=95%,Q30>=85%
比对至参考基因组
软件:Hisat2/STAR
Input:Clean FASTQ
Output:SAM
选择spliced aligner
Hisat2(http://daehwankimlab.github.io/hisat2/)
(https://www.jianshu.com/p/ce3f4afb9b60)
(https://www.jianshu.com/p/479c7b576e6f)
STAR(https://www.jianshu.com/p/5b6dfc954315)
比对文件过滤
软件:samtools(去除低质量比对)/sambamba/picard
Input:Raw SAM
output:Clean BAM
去除低质量比对(MAPQ<30,MAPQ)
去除多重比对(一条read比对到基因组的多个地方)
去除duplication(不同的reads比对到基因组的同一位置)
去除比对至black list region (ENCODE数据库,大量reads比对到的区域)的reads
MAPQ可信度:mapping的质量值,错误率的-10log10(https://www.jianshu.com/p/9c87bba244d8 ,https://zhuanlan.zhihu.com/p/35495052)
samtools(http://samtools.sourceforge.net/)
(https://www.jianshu.com/p/6b7a442d293f)
sambamba(http://lomereiter.github.io/sambamba/)
picard(http://broadinstitute.github.io/picard/)
RNA甲基化生信分析
鉴定m6A修饰
输入bam文件后就可以拿到peak
软件:exomePeak/MeTPeak
Input:BAM
Output:peak
exomePeak(https://bioconductor.riken.jp/packages/3.1/bioc/html/exomePeak.html)
MeTPeak(https://github.com/compgenomics/MeTPeak)
鉴定差异m6A修饰
差异基因有峰(IGV可视)
软件:exomePeak(最好有重复的样本)/MeDiff(?)/QNB(有重复或不重复数据)/RADAR
Input:BAM
Output:差异peak
QNB(https://rdrr.io/cran/QNB/)
RADAR(https://scottzijiezhang.github.io/RADARmanual/workflow.html)
Peak/差异peak注释及可视化
软件:ChIPpeakAnno/ChIPseeker
Input:BAM
OutPut:差异peak
样本间peak的重叠情况(ChIPseeker)
Peak在基因元件上的分布
Peak在基因元件上的分布模式(Guitar)
特定基因/区域的peak(BAM->TDF->IGV)
注意:最好注释到基因水平
m6A-seq建库用的片段小,无法注释到转录本水平,最好注释到基因水平
需区分正负链
ChIPpeakAnno(https://www.jianshu.com/p/9b89c12f79d0)
ChIPseeker(https://www.jianshu.com/p/c76e83e6fa57)
Peak/差异peak关联基因的功能富集
软件:AllEnricher/clusterProfiler/ReactomePA
Input:GeneSymbol/Entrez ID,相应数据库
output:富集结果及图片展示
AllEnricher(https://github.com/zd105/AllEnricher)
clusterProfiler(https://www.jianshu.com/p/d484003dced5)
ReactomePA(https://cloud.tencent.com/developer/article/1695223)
Peak和差异peak的motif鉴定
软件:HOMER/MEME(可在线分析)
Input:peak/sequence(FASTA)
output:de novo motif、known motif
注意:区分正负链(-rna参数,根据正负链进行分析)
HOMER(https://zhuanlan.zhihu.com/p/73504094)
MEME(https://jingyan.baidu.com/article/154b46318a419b69ca8f41c5.html)
与RNA-seq关联分析
Input+IP:m6A修饰图谱;
Input:表达图谱 Input文库可当RNA-seq
四象限图,m6A、RNA差异表达。
学徒作业
学习He lab m6A-seq_analysis_workflow,链接是:https://github.com/scottzijiezhang/m6A-seq_analysis_workflow
无非就是hisat2的比对,拿到bam文件后就可以直接call peaks,拿到bed格式的peaks坐标文件后就各种下游分析,基因注释或者motif查找,以及部分可视化:
比较特殊的是:R package MeRIPtools for peak calling. homer软件的findMotifs.pl 程序,针对bedtools getfasta拿到的peaks坐标附件序列文件 可视化,差异分析,功能注释,类比于其它表观组学
自己摸索好这个流程,然后完成发表在2020年11月的文章:LncRNA DILA1 inhibits Cyclin D1 degradation and contributes to tamoxifen resistance in breast cancer 的实战。